;; Various series, most prominently Taylor and Maclaurin.

(in-package :azala)

(defun taylor (exp dx a degree)
  "Return the Taylor series approximation for an expression to a given
degree, centered on a, with respect to dx."
  (variable-arity-to-binary
   '+
   (loop for n from 0 to degree
	 collect `(/ (* ,(substitute-variable (diff-degree exp dx n)
					      dx a)
		      (expt (- ,dx ,a) ,n))
		   ,(! n)))))

(defun variable-arity-to-binary (operator arguments)
  (simplify (reduce #'(lambda (x y)
			(list operator x y))
		    arguments :from-end t)))

(defun maclaurin (exp dx degree)
  "Find a Maclaurin series approximation. This is just TAYLOR with one
parameter omitted."
  (taylor exp dx 0 degree))

(declaim (inline maclaurin))

(maclaurin '(sin x) 'x 9)
(maclaurin '(cos x) 'x 9)

(maclaurin '(expt %e x) 'x 9)

(defmacro def-taylor (name args exp dx a degree)
  `(defun ,name ,args
    (declare (double-float ,@args))
    (the double-float
      ,(taylor exp dx a degree))))

(def-taylor my-cos (x)
  (cos x) x 0 9)

(def-taylor my-sin (x)
  (sin x) x 0 9)

(def-taylor my-e^x (x)
  (expt %e x) x 0 15)

(def-taylor my-e^2x (x)
  (expt %e (* x 2)) x 0 15)

;(def-taylor badfunc (x)
;  (/ 1 (- 1 (expt x 7))) x 0 10)
;
;(maclaurin '(/ 1 (- 1 (expt x 7))) 'x 10)